library(NMF)


try({
  getwd()
  setwd('./一致性聚类/')
})

BiocManager::install()

gene_exp <- read.table("m7Gexp.txt")
gene_exp <- gene_exp[-1,]

gene_exp <- as.matrix(gene_exp)
gene_exp <- apply(gene_exp, 2, as.numeric)
x <- log(gene_exp + 1) 
x.rank <- nmf(x,2:10,nrun=10,seed=111111,method = 'brunet') 
# 参数解释：file，rank=2:10（rank在2-10之间），nrun=10（迭代次数为10次），method选择分析方法
summary(x.rank)
pdf('nmfxrank.pdf',width = 15,height = 15,onefile = F)
plot(x.rank) #选择最佳rank值，最大值或波动最大处的最高点
dev.off()
#使用最佳rank，或者最符合需求的rank
x.rank2 <- nmf(x,5,nrun=10,seed=555555,method = 'brunet') # run parameters,rank=5
# coefmap from multiple run fit: includes a consensus track


##混合系数矩阵：coefmap 作用：展示从最佳拟合结果中获得的簇（聚类数）和一致性矩阵的层次聚类。
pdf('coefmap_rank3.pdf',width = 15,height = 15,onefile = F)
par(cex = 3)
coefmap(x.rank2,
        annRow = NA,
        annCol = NA,
        main = "Metagene contributions in each sample",
        info = FALSE)
dev.off()


##基底矩阵：basismap 作用：每一行显示主导的基底组分，即每一行有最高负载的基底组分。
pdf('basismap_rank3.pdf',width = 15,height = 15,onefile = F)
basismap(x.rank2,
         annRow = NA,
         annCol = NA,
         main = "Metagenes",
         info = FALSE)
dev.off() 
##一致性矩阵：consensusmap 作用：基于指定rank评估聚类稳定性的方法，考虑由多个独立NMF运行结果计算得到的连接矩阵。
pdf('consensusmap_rank3.pdf',width = 15,height = 15,onefile = F)
consensusmap(x.rank2,
             annRow = NA,
             annCol = NA,
             main = "Consensus matrix",
             info = FALSE)
dev.off() 
##查看分组情况
# group <- predict(x.rank2)
# group <- as.data.frame(group)
# group$group <- paste0('Cluster',group$group)
# group$sample <- rownames(group)
# group<- group[order(group$group),]
# table(group$group)
# head(group)
# #保存分组情况
# save(group,file = 'expr_group3.rdata')  #建立稳定模型后，可以直接打开rdata对数据进行分类分析
# write.csv(group,'expr_group3.csv') #也可以以csv格式导出到本地，查看具体分类情况
